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We will describe how a new, quite simple, but highly effective algorithm, together with the asymptotically fast FFT-based 
high-precision number multiplication of Mathematica 4 can calculate the ground state of the x 4 anharmonic oscillator to the 
new record of more than 1000 digits. 



■ Introduction 



As it is well known, only a very limited number of one-dimensional potentials allow for an exact solution of the 
Schrodinger equation. This means that for many model potentials we must resort to numerical solution methods. For 
judging their accuracy, reliability, and speed, it is important to have high-precision values of certain nonexactly solvable 
potentials. The most investigated of such potentials is the quartic anharmonic oscillator (see [1] to [19]), described by 



-^'(z) + z 4 <Mz) = £*^(z) CD 
The eigenfunctions to the eigenvalues e k decay exponentially for z -> ± oo . 



■ The Hill determinant method 



A classical method for solving Sturm-Liouville problems of type (1) is to calculate the eigenvalues of a truncated version of 
the corresponding Hill determinant. Using the harmonic oscillator basis <f>„(z) we write tffoiz) = Ti^Lo a k Qkiz) where 

-<P'liz) + z 2 <f> k (z) = e k <t> k (z) (2) 



e~- H n (z) (3) 



V V?r 2" n ! 

Thus we form the matrix elements h mn = f <p m (z) (-((>„ (z) + z 4 4>n (z)) dz . For n > m we obtain 



n - m > 4 



2T i - 4 y (32 n 5 m ,„- 2 (n - l) 2 + 16 (n - 3) (n - 2) n 5 m ,„_ 4 (n - 1) + else (5) 
4 (2 n (3 n + 5) + 5) S n , n + 8 (n + 1) 6 m ,„ +2 + 6 m , n+ 4 ) 



A rough estimation shows that we obtain about 0.2 digits per harmonic oscillator state. So, by taking into account the first 
500 eigenstates and carrying out the calculation with about five thousand digits we obtain about 120 reliable digits for en . 
(This calculation takes about 20 minutes on a year-2000 vintage workstation using Mathematica 4 [20].) 



£ = 1 .0603620904841 8289964704601 66926635455 1 52087285289779332 1 624524 1695 
943563044344421126896299134671703510546244358582525580982763829... 

The Hill determinant approach allows the calculation of the eigenvectors in addition to the calculation of the eigenvalues. 
The following graphic visualizes the matrix of eigenvectors of {h mA ) x< <1Q0 . The graphic shows that the lowest eigenfunc- 
tions are quite similar to the harmonic oscillator eigenfunctions. Higher states are complicated mixtures of harmonic 
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oscillator states. The overall "checkerboard"-like structure results from the fact that the contribution of the antisymmetric 
(symmetric) harmonic oscillator states to the symmetric (antisymmetric) anharmonic oscillator states is identically zero. The 
very high states are dominated by truncation effects and do not correctly mimic the anharmonic oscillator states. 
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■ The new algorithm 

To get a very high-precision approximation of 

-i//'(x) + zV« = <^(x) (6) 
we start with the series expansion 

/; 

^(x)=yn(x) = Y i ^W^ (7) 
/<=() 

For the ground state we choose (ignoring normalization) i/f(0) =1, ^'(0) = 0. For "suitable chosen" x* we then find 
high-precision approximations for the zeros of y„ (x* ) and y' n (x* ) . These zeros then bound Ao from below and above. 

Using the differential equation, one obtains the following recursion relation for (A) : 

Om-6(A) - Aa m _2(A) 

a* (A) = ; (8) 

mr — m 

For large n (n -» oo ) we want the function y n (x) to vanish as x -» oo . For a A smaller than the exact eigenvalue, the function 
y n (x* ) will not have a zero, but the function y' n (x* ) will have a zero for a certain x* . For a A larger than the exact eigen- 
value, the function y n (x* ) will have a zero, but the function y' n (x* ) will not have a zero for a certain x* . This fact allows us 
to find a bounding interval for Ao . The next two graphics show yso(x) and y'go(x) for 10 equidistant values for A from the 
interval [1.05, 1.08] to visualize this bounding process. (For more details, see [22].) 
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It is straightforward to implement the calculation of the bounding interval for Ao in Mathematica in a three-line program 
(see [21]). Using FindRoot we calculate high-precision values for the zeros of y„ (£) and y' n (£) . 

ABounds[n_, f_, opts ] : = Map [Function [f, 

A /. FindRoot [f[n, X, f] ==0, {A, 106/ 100, 107 / 100}, opts]], {y, yPrime}]; 
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The calculation of y„ (£) and y' n (£) is also straightforward, based on recursive calculations of a k (A) . 

y[n_, A_Real, £_] : = Module [ {a6, a4, a2, ak, ct} , 

{a6, a4, a2} = {1, -A / 2, A 2 / 24}; a = a6 + a4 f 2 + a2 f 4 ; 
Do[ak = a6 - Aa2 / (k (k - 1) ) ; {a6, a4, a2} = {a4, a2, ak} ; 
a= a + ak£ k , {k, 6, n, 2}]; a] 

y[n_, A_Real, f_] : = Module [ {a6, a4, a2, ak, a), 

{a6, a4, a2} = {1, -A/ 2, A 2 / 24}; ct = a6 + 2 f a4 + 4 f 3 a2 ; 
Do[ak = a6 - Aa2 / (k (k - 1) ) ; {a6, a4, a2} = {a4, a2, ak} ; 
a= otkakf 1 " 1 , {k, 6, n, 2}]; a] 

Now we calculate ABounds [ 1 6000 , 16, startingValues , WorkingPrecision^6000, AccuracyGoal^600 , 
Maxiterations^lOO] (where startingValues has been obtained from a call to ABounds recursively). In a few minutes, we 
get a 1 184-digit approximation to the ground state energy of the quartic anharmonic oscillator. 

£o = 

1 .060362090484 1 8289964704601 66926635455 1 52087285289779332 1624524 1 69594356304434442 1 1 26896299 1 34671703'-. 
51054624435858252558087980821029314701317683637382493578922624600470817544696014163748841728225690: 
93575779088806178879026360154939569027519614890094293487358440944269489790121397146429095192335453: 
8283470335057576151120257039888523720240221841103086573731091398915453658410311167940583354860009Z 
744006963112670238862297142969961059215583226671376935508673610000831830027517926233573913906136181 
77649859696181499412792809272840707956106044072294680994913627572927387279136890279842472226171694 
48895475137043806840543918778772953234245874372543178323190603810687416044034374530146847278139186 
294047043103401351071607110353008929823275427661518986950565047160252756089526262191025688200964411 
28781564005270529293240507638265028259112477362538471854714402572285438485297450458570978840249066' 
99570476844587709 1 762029 124375273254907 1 16433440230294730692398 1908956853745359884460 1 60023 1 329 193: 
05939586930491664428163394616332428700426146123774300995223420420859773569015356541685030894185134! 
795734106585479719467596466796613467688586437952654519560568286715958338884743467012042420714919291 
048732... 

Statistical analysis of the number does not show any regularity. 



■ Summary 

A power series-based approach to the high-precision calculation of the ground state of the anharmonic oscillator was 
presented. Mathematica code to carry out the calculation, as well as results, were given. The method can straightforwardly 
be used to calculate tens of thousands of digits of the quartic anharmonic, as well as other anharmonic oscillators. Work 
concerning the application of the method to higher states is in progress. 

All calculations and visualizations have been carried out in Mathematica 4. 
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